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Abstract 

We  perform  calculations  of  2D  laminar  and  3D  turbulent  channel  flow  with 
periodic  streamwise  boundary  conditions.  Prom  the  laminar  results  we  verify 
that  the  Fire  Dynamics  Simulator  (FDS)  wall  boundary  condition  is  second-order 
accurate.  For  the  turbulent  cases  we  adapt  the  wall  stress  model  of  Werner  and 
Wengle  to  FDS  and  verify  the  implementation  by  replicating  the  Moody  diagram, 
a  plot  of  friction  factor  versus  Reynolds  number  for  turbulent  pipe  flow. 


1     Introduction 

Wall  flows  are  notoriously  challenging  for  large-eddy  simulation  (LES)  [3,  4,  5,  10,  11].  In 
spite  of  their  promise  and  sophistication,  practical  LES  codes  are  resigned  to  model  the  wall 
shear  stress  as  opposed  to  resolving  the  dynamically  important  length  scales  near  the  wall. 
In  this  work  we  introduce  the  Werner  and  Wengle  (WW)  wall  model  [13]  into  the  National 
Institute  of  Standards  and  Technology  (NIST)  Fire  Dynamics  Simulator  (FDS)  as  a  practical 
first  step  in  developing  models  for  turbulent  flow  around  complex  geometry  and  over  complex 
terrain.  Such  models  are  required  in  order  for  FDS  to  accurately  model,  for  example,  tunnel 
fires,  smoke  transport  in  complex  architectures,  and  wildland-urban  interface  (WUI)  fires 
[1].  As  a  minimum  requirement,  a  wall  model  should  accurately  reproduce  the  mean  wall 
stress  for  fiow  in  a  straight  channel.  We  verify  that  this  is  true  for  FDS  by  reproducing  the 
Moody  chart,  a  plot  of  friction  factor  versus  Reynolds  number  for  pipe  flow  [8]. 


*  Email:  raiidall.incdermott@nist.gov 
NIST  Technical  Note  1640  July  2,  2009 


The  remainder  of  this  article  is  organized  as  follows.  In  Section  2  we  describe  the  model 
formulation.  In  Section  3  we  give  an  overview  of  the  WW  model.  Then,  in  Section  4,  we 
conduct  a  verification  study  of  the  wall  boundary  conditions  for  laminar  and  turbulent  flows 
in  FDS.  From  this  study  we  are  able  to  draw  quantitative  conclusions  in  Section  5  about  the 
accuracy  of  the  channel  flow  simulations.  A  derivation  of  the  FDS  implementation  of  the 
WW  model  is  given  in  Appendix  A  and  the  friction  factor  for  2D  Poiseuille  flow  is  derived 
in  Appendix  B. 

2     Formulation 

Details  of  the  FDS  formulation  are  given  in  the  Technical  Guide  [7] .  Here  we  provide  only 
the  salient  components  of  the  model  necessary  for  treatment  of  constant  density  channel 
flow. 

The  filtered  continuity  and  momentum  equations  are: 


dui      duiUj  1 


-^  +  ^  +  ?Zii  +  ^^L- 


dt         dxj  p  [dxi      dxi      dxj        dx 


(2) 


where  r*,?*  =  p{uiU]  —  UiUj)  is  the  subgrid-scale  (sgs)  stress  tensor,  here  modeled  by  gra- 
dient diffusion  with  dynamic  Smagorinsky  [6]  used  for  the  eddy  viscosity.  In  this  work  we 
specify  a  constant  pressure  drop  dp/dx  in  the  streamwise  direction  to  drive  the  flow.  The 
hydrodynamic  pressure  p  is  obtained  from  a  Poisson  equation  which  enforces  (1). 

When  (2)  is  integrated  over  a  cell  adjacent  to  the  wall  in  an  LES  it  turns  out  that  the 
most  difficult  term  to  handle  is  the  viscous  stress  at  the  wall,  e.g.  fxz\z=o,  because  the  wall- 
normal  gradient  of  the  streamwise  velocity  component  cannot  be  resolved.  Note  that  the  sgs 
stress  at  the  wall  is  identically  zero.  We  have,  therefore,  an  entirely  different  situation  than 
exists  in  the  bulk  flow  at  high  Reynolds  number  where  the  viscous  terms  are  negligible  and 
the  sgs  stress  is  of  critical  importance.  The  quality  of  the  sgs  model  still  influences  the  wall 
stress,  however,  since  other  components  of  the  sgs  tensor  affect  the  value  of  the  near- wall 


velocity  and  hence  the  resulting  viscous  stress  determined  by  the  wall  model.  In  particular, 
it  is  important  that  the  sgs  model  is  convergent  (in  the  sense  that  the  LES  formulation 
reduces  to  a  DNS  as  the  filter  width  becomes  small)  so  that  as  the  grid  is  refined  we  can 
expect  more  accurate  results  from  the  simulation. 

The  model  used  for  t^j  =  fxz\z=o  in  this  work  is  the  Werner  and  Wengle  model  [13]  which 
we  describe  in  more  detail  below. 

3     The  Werner  and  Wengle  Wall  Model 

An  important  scaling  quantity  in  the  near-wall  region  is  the  friction  velocity,  defined  as 


u*  =  y/r^/p.  Prom  the  friction  velocity  we  define  the  nondimensional  streamwise  velocity 
u'^  =  u/u*  and  nondimensional  wall- normal  distance  z'^  =  z/i,  where  £  =  ii/{pu*).  The  law 
of  the  wall  is  then  given  by  [10,  12] 

u+    =    2+  for     2+  <  5  (3) 

w+    =    2.4  In  2+ +  5.2    for     2+>30,  (4) 

The  region  5  <  ^^  <  30,  where  both  viscous  and  inertial  stresses  are  important,  is  referred 
to  as  the  buffer  layer.  The  upper  range  of  the  log  law  depends  on  the  Reynolds  number 
[10,  14]. 

Werner  and  Wengle  [13]  propose  a  simplification  to  the  law  of  the  wall  to  eliminate  the 
mathematical  difficulties  of  handling  the  buffer  and  log  layers.  Furthermore,  WW  suppose 
that  their  simplified  formula  for  the  streamwise  velocity  holds  instantaneously  within  the 
LES.  The  WW  wall  profile  is  [11] 

u+    =    2+  for     2+  <  11.81  (5) 

u+    =    A{z+f    for     z+  >  11.81,  (6) 

where  ^  ==^  8.3  and  B  =  1/7.  Note  that  a  power  law  has  been  substituted  for  the  log  law 
and  the  viscous  sublayer  and  the  power  law  region  are  matched  within  the  buffer  region.  A 
comparison  of  the  log  law  and  the  power  law  is  shown  in  Figure  1.   In  the  region  11.81  < 
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Figure  1:  The  law  of  the  wall.  We  have  omitted  the  buffer  layer  since  it  is  not  considered  in  the 
WW  model.  For  2"*"  ^  11.81  we  have  the  viscous  sublayer.  For  z^  >  11.81  we  show  a  comparison 
of  the  log  law  (4)  (dashed  line)  and  the  WW  power  law  (6)  (solid  line)  with  A  =  8.3  and  B  =  1/7. 


z'^  <  10^  the  power  law  is  a  good  approximation  to  the  log  law  and  for  2"*"  >  10^  the  power 
law  loosely  exhibits  wake  region  behavior  for  a  flow  with  Re  w  5  x  10^  [10,  14].  As  we  see 
below,  this  functional  behavior  has  consequences  for  high  Re  flows. 

For  the  purposes  of  adapting  the  WW  model  to  FDS  we  suppose  that  the  first  off-wall 
velocity  component  u  represents  the  WW  profile  averaged  in  the  wall-normal  direction  (refer 
to  Figure  2).  The  density  is  taken  as  the  average  of  the  neighboring  cell  values  and  uniform 
along  the  face.  The  WW  model  as  implemented  in  FDS  is  then  given  by 
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Figure  2:  Near- wall  grid. 
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The  parameters  defined  in  (9)-(12)  are  introduced  to  provide  a  one-to-one  correspondence 
between  this  paper  and  the  FDS  code.  Note  that  p,  is  the  average  of  the  molecular  viscosity 
from  the  neighboring  cells.  A  detailed  derivation  of  (8)  is  given  in  Appendix  A. 

In  order  to  decide  which  formula  to  use  for  the  wall  stress,  (7)  or  (8),  we  must  know  z'^, 
which  of  course  depends  on  t^j.  As  a  practical  matter  of  implementation,  given  that  most 
boundary  layers  in  FDS  are  under-resolved,  we  first  calculate  r^,  from  (8);  we  then  obtain 
z^  —  \ItwIp  ;  if  -z"*"  >  11.81,  then  the  computed  value  of  r^,  is  retained,  else  r^^  is  taken 
from  (7),  which  actually  involves  no  additional  computation  since  the  ghost  cell  value  for 
the  velocity  is  prescribed  for  a  no-slip  wall  by  default. 


4     Results 

4.1  Laminar 

As  verification  of  the  no-slip  boundary  condition  and  further  verification  of  the  momentum 
solver  in  FDS,  we  perform  a  simple  2D  laminar  (Poiseuille)  flow  calculation  of  flow  through 
a  straight  channel.  The  FDS  input  files  are  stored  in  the  repository  [2]  under  poiseuille_*. 
The  height  of  the  channel  is  i:r  =  1  m  and  the  length  of  the  channel  is  L  =  8  m.  The  number 
of  grid  cells  in  the  streamwise  direction  x  is  N^  —  8.  The  number  of  cells  in  the  wall-normal 
direction  z  is  varied  Nz  =  {8,16,32,64}.  The  fluid  density  is  p  =  1.2  kg  m""^  and  the 
viscosity  is  0.025  kg  m~^  s~^.  The  mean  pressure  drop  is  prescribed  to  be  dp/dx  =  —  1  Pa 
m"^  resulting  in  Re/f  «  160.  The  (Moody)  friction  factor  /,  which  satisfies 

Ap  =  f^^-pu\  (13) 

is  determined  from  the  steady  state  mean  velocity  u  which  is  output  by  FDS  for  the  specified 
pressure  drop.  The  exact  friction  factor  for  this  flow  is  /exact  =  24/Re//  (see  Appendix  B). 
The  firiction  factor  error  \f  —  fexactl  is  plotted  for  a  range  of  grid  spacings  Az  =  H/Nz  in 
Figure  3  demonstrating  second-order  convergence  of  the  laminar  velocity  field. 

4.2  Turbulent 

To  verify  the  WW  wall  model  for  turbulent  flow  we  perform  3D  LES  of  a  square  channel  with 
periodic  boundaries  in  the  streamwise  direction  and  a  constant  and  uniform  mean  pressure 
gradient  driving  the  flow.  The  calculation  is  nearly  identical  to  the  laminar  cases  of  the 
previous  section,  except  here  we  perform  3D  calculations  and  maintain  cubic  cells  as  we 
refine  the  grid:  we  hold  the  ratio  8:1:1  between  Nx.Ny-.Nz  for  all  cases.  The  cases  shown 
below  are  identified  by  their  grid  resolution  in  the  z  direction.  The  velocity  field  is  initially 
at  rest  and  develops  in  time  to  a  mean  steady  state  driven  by  the  specified  mean  pressure 
gradient.  The  presence  of  a  steady  state  is  the  result  of  a  balance  between  the  streamwise 
pressure  drop  and  the  integrated  wall  stress  from  the  WW  model.  FDS  outputs  the  planar 
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Figure  3:  FDS  exhibits  second-order  convergence  for  laminar  (Poiseuille)  flow  in  a  2D  channel. 


average  velocity  in  the  streamwise  direction  and  once  a  steady  state  is  reached  this  value 
is  used  to  compute  the  Reynolds  number  and  the  friction  factor.  Table  1  provides  a  case 
matrix:  nine  cases  for  three  values  of  specified  pressure  drop  and  three  grid  resolutions.  The 
nominal  Reynolds  number  (obtained  post-run)  is  listed  along  with  the  friction  factor  from 
the  most  refined  FDS  case  and  the  friction  factor  computed  (iteratively)  from  the  Colebrook 
equation, 
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77 


-  -2.0  log 


10 


e/D        2.51  \ 
rT  ^  Rev?/ 


(14) 


which  is  a  fit  to  the  turbulent  range  of  the  Moody  chart  (see  e.g.  [9]).  The  parameter  s/D 
is  the  relative  roughness  where  D  is  the  hydraulic  diameter  of  the  pipe  or  channel  and  Re  is 
the  Reynolds  number  based  on  D.  For  all  the  cases  reported  here  the  hydraulic  diameter  is 
equivalent  to  the  channel  height,  D  =  H,  and  the  walls  of  the  channel  are  smooth,  i.e.  e  =  0. 
FDS  input  files  are  stored  in  the  repository  [2]  as  moody_*.  To  provide  a  qualitative  picture 
of  the  flow  field,  Figure  4  shows  contours  of  streamwise  velocity  for  the  case  dp/dx  =  —  1  Pa 


Table  1:  Case  matrix  and  friction  factor  results  for  turbulent  channel  flow.  The  height  of  the  first 
grid  cell  Az  is  given  in  viscous  units  z'^  for  each  case.  Additionally,  the  table  gives  the  nominal 
Reynolds  number  Re//  and  the  FDS  friction  factor  results  compared  to  the  Colebrook  equation 

(14). 


dp/dx 

2+ 
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/FDS 

/  Colebrook 

rel.  error 

(Pa/m) 
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0.0081 

6.0 

m-i  and  N^  =  32. 

In  Figure  5  we  replicate  the  Moody  diagram  for  both  the  laminar  and  turbulent  cases 
presented  in  this  work.  The  laminar  cases  provide  the  exact  result  for  two  different  Reynolds 
numbers  and  the  turbulent  cases  are  converging  to  the  empirical  values  of  the  friction  factor 
for  smooth  pipes.  It  is  interesting  to  compare  the  turbulent  results  in  Figure  5  with  the  values 
of  2+  shown  in  Table  1.  Notice  from  the  table  that  the  z"*"  values  for  the  dp/dx  =  —100  cases 
are  in  a  range  where  the  power  law  (6)  deviates  significantly  from  the  log  law  (4)  (see  Figure 
1)  and  this  may  explain  why  the  results  for  this  high  Re  cgise  are  somewhat  grid  sensitive. 

5     Conclusions 

In  this  work  we  have  verified  the  FDS  wall  model  for  both  laminar  and  turbulent  flow 
through  straight  channels.  We  have  shown  that  for  the  laminar  (DNS)  case  FDS  is  second- 
order  accurate.  It  is  suggested  elsewhere  that,  as  a  rule  of  thumb,  10  %  accuracy  is  the 
best  that  can  be  expected  from  friction  factor  calculations  of  turbulent  flow  [9].  We  have 
adapted  the  Werner  and  Wangle  wall  model  to  variable  density  flows  (though  only  constant 
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Figure  4:  LES  of  square  channel  flow  with  smooth  walls  and  periodic  streamwise  boundaries  using 
dynamic  Smagorinsky  and  the  Werner  Wengle  wall  model.  For  this  image  N^  =  32  and  the  mean 
pressure  drop  is  dp/dx  =  —  1  Pa  m~^  resulting  in  Ren  =  7.5  x  10^  and  a  friction  factor  of  /  =  0.0128. 
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Figure  5:  The  FDS  Moody  Chart.  The  sohd  hne  for  Re  <  2000  is  the  analytical  result  for  2D 
Poiseuille  flow,  /  —  24/Re  (see  Appendix  B).  The  sohd  lines  for  Re  >  2000  (from  the  Colebrook 
equation  (14))  are  for  turbulent  flow  at  various  levels  of  relative  roughness  e/D  shown  on  the  right 
axis.  Stars  are  DNS  results  from  FDS  at  a  single  grid  resolution  (Nz  —  64)  and  the  symbols  are 
FDS  results  for  3D  LES  with  dynamic  Smagorinsky  and  the  Werner  Wengle  wall  model  at  three 
grid  resolutions  {N^  =  {8, 16, 32}). 
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density  flows  are  tested  here)  for  smooth  walls  and  have  shown  that,  when  combined  with 
the  dynamic  model  for  the  eddy  viscosity,  FDS  is  capable  of  reproducing  friction  factors  for 
a  broad  range  of  Reynolds  numbers  to  within  6.0  %  relative  accuracy. 
In  future  work  we  plan  to: 

(1)  Incorporate  roughness  effects  into  the  wall  model  for  both  constant  density  and  stratified 
flows  relevant  to  atmospheric  turbulence. 

(2)  Increase  the  accuracy  of  the  FDS  immersed  boundary  method  for  flow  around  complex 
geometry  (i.e.  eliminate  the  requirement  of  "stair-stepped"  geometry).  In  particular,  the 
results  for  the  flows  presented  here  should  be  invariant  to  a  rotation  of  the  computational 
mesh  relative  to  the  streamwise  direction  of  the  flow.  At  present,  this  is  not  the  case  in 
FDS. 

(3)  Validate  the  FDS  convective  heat  transfer  models  against  Nusselt  number  correlations 
for  Rayleigh-Benard  convection. 

A     Derivation  of  the  WW  Model 

To  obtain  (8)  we  take  the  first  off-wall  value  of  the  streamwise  velocity  to  be 
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and  then  substitute  the  WW  profile  for  u{z)  and  integrate. 

Let  Zm  denote  the  dimensional  distance  from  wall  where  z'^  =  11.81.    Equation  (15) 
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We  will  integrate  /  and  //  separately.  First,  however,  we  must  find  a  way  to  eliminate 
the  unknown  Zm-  To  do  this  we  equate  (5)  and  (6)  at  the  point  where  the  viscous  and  power 
law  regions  intersect,  i.e.  z'^  =  11.81  =  2+  =  z-mpu* /jl. 
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We  now  have  Zjn  in  terms  of  r^  and  otherwise  known  values. 
Integrating  section  7  of  (16)  we  find 
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Integrating  section  //  yields 
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Plugging  (18)  and  (19)  back  into  (16)  gives 
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Rearranging  for  t^  we  find 
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'^%i±i?f^,  ^ 


Tw      ^      P 


1    —   5     ,i+B     /        /i 


2  \pAz 

which  corresponds  to  Eq.  (9.46)  in  [11]. 


l+B 


+ 


A      \pAz 
l  +  B  f    p 


A      \pAz 


B     1 

u 


l  +  B 


(19) 


(20) 


(21) 


B     Poiseuille  Flow  in  2D 


We  consider  fully  developed  laminar  flow  in  a  2D  straight  channel  of  length  L  and  height  H. 
The  momentum  balance  in  the  streamwise  direction  gives 


dp         d^u 
dx         dz"^ 
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(22) 


Integrating  (22)  over  the  height  of  the  channel  we  find 


^^^^-ifj^'-"^^-  (23) 


Prom  (23)  we  find  the  average  streamwise  velocity  is 

1    f^ 
u  =  —        u{z)  dz 

2/i  dxJo^  ^ 


'^"W  (24) 


1     ^Ptt2 


Vljidx 
The  (Moody)  friction  factor  /  is  defined  such  that  the  pressure  drop  Ap  is  given  by  [9] 

/^P  =  f^\pu\  (25) 


Combining  (24)  with  (25)  and  noting  that  ^  =  -^  we  find 


A.  =  /|i. 


1    Ap^, 


12^  L 


u,  (26) 


and  finally 


24^^J^  (27) 

•'       pHu      Ren'  ^     ^ 
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